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We investigate the low-lying eigenmodes of the Dirac matrix with the aim to gain more insight 
into the temperature dependence of the anomalous Ua{ 1) symmetry. We use the overlap operator 
to probe dynamical QCD configurations generated with (2+l)-flavors of highly improved staggered 
quarks. We find no evidence of a gap opening up in the infrared region of the eigenvalue spectrum 
even at 1.5 T c , T c being the chiral crossover temperature. Instead, we observe an accumulation of 
near-zero eigenmodes. We argue that these near-zero eigenmodes are primarily responsible for the 
anomalous breaking of the axial symmetry still being effective. At 1.5 T C1 these near-zero eigenmodes 
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temperature, the average size of the instantons is 0.223(8) fm and their density is 0.147(7) frn -4 . 
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I. INTRODUCTION 


Owing to the near-degeneracy and smallness of the 
up and down quark masses, the Quantum Chromo¬ 
dynamics (QCD) Lagrangian possesses an approximate 
U l {2) x £7^(2) = SU L { 2) x SU R ( 2) x U v ( 1) x U A ( 1) 
symmetry. The fact that we do not see parity doublet 
hadrons in our world implies that the SUl{ 2) x SU R ( 2) 
chiral symmetry is spontaneously broken down to the 
SUy (2) isospin symmetry of the vacuum. It is well known 
from first principle lattice QCD studies that above the 
chiral crossover M temperature of T c = 154(9) MeV 
0 the chiral symmetry of QCD gets restored. 

On the other hand, the axial U A (1) symmetry of the 
QCD Lagrangian is always broken due to the presence of 
quantum fluctuations. This gives rise to the well-known 
anomalous non-conservation of the axial current 10. 
The explicit violation of the global U A { 1) symmetry is 
due to the presence of topologically nontrivial gauge field 
configurations [§}. Although U A ( 1) is not an exact sym¬ 
metry of QCD, the magnitude of its breaking near T c is 
expected to influence the nature of the chiral phase tran¬ 
sition in the limit of two vanishingly small light quark 
masses. Perturbative renormalization group studies of 
model quantum field theories with the same global sym¬ 
metries as QCD suggest that if U A (l) is not effectively 
restored at T c , the chiral phase transition is of second 
order, belonging to the 3-dimensional 0(4) universality 
class (Mi- If the axial symmetry gets effectively re¬ 
stored for T ~ T ci _the chiral phase transition can be ei¬ 
ther of first order 0, [HI or of second order with the sym¬ 
metry breaking pattern Ul{ 2) x Ur(2) —> Z7y(2) fillll^| . 
In order to resolve the nature of the phase transition of 
QCD with two light quark flavors it is thus important to 
understand the significance of the anomalous U A ( 1) in 
the high temperature phase. 

At low temperatures, 17,4(1) is also broken explicitly 
by the presence of a non vanishing vacuum chiral con¬ 
densate. In the chirally symmetric phase, the vacuum 


condensate vanishes and the mechanism of global 17,4(1) 
breaking can be studied directly. The microscopic mech¬ 
anism for U A (1) breaking in the chirally symmetric phase 
of QCD presents an intriguing puzzle. The chiral conden¬ 
sate, which is the order parameter related to the restora¬ 
tion of chiral symmetry in QCD with massless quarks, 
can be expressed in terms of the eigenvalues A of the 
Dirac operator as 
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where p( A, m) is the eigenvalue density. On the other 
hand, U A (1) is not a global symmetry, so one cannot 
define a corresponding order parameter. For two light 
quark flavors, an approximate restoration of U A (1) would 
result in the degeneracy of the correlation functions of the 
pion and the scalar iso-triplet delta meson [Q. Specifi¬ 
cally, the difference of the integrated correlation functions 
of these mesons in terms of the eigenvalues of the Dirac 
operator is 
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In the limit of vanishingly small quark mass m and infi¬ 
nite volume V, the chiral condensate is proportional to 
the density of near-zero eigenvalues in accordance with 
the Banks-Casher relation & (ifrtp) = 7rp(0,0). For 
T >T C , chiral symmetry gets restored and the chiral or¬ 
der parameter (i/np) vanishes implying that p{ 0, 0) must 
also vanish. Motivated by the free theory limit at fi¬ 
nite temperature, where the spectral density p( A, 0) has 
a gap up to the lowest fermion Matsubara frequency, 
0 < A < 7rT, one possibility by which chiral symmetry 
restoration may occur in the chiral limit is through the 
generation of a gap in the infrared part of the eigenvalue 
spectrum. Such a scenario, however, would also lead to 
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the vanishing of x-k ~ Xs, he. to the effective restoration 
of both chiral and Ua(1) symmetry. 

In fact, more rigorous calculations based on chiral 
Ward identities for up to 4-point correlation functions 
show []1| that if the eigenvalue density for QCD with 
two light quark flavors is an analytic function in m 2 , 
it must have the form lim m _>o /o(A, to) ~ A 3 + 0( A 4 ) in 
the chirally symmetric phase, similar to that for the free 
theory at T = 0. It was further shown [TEJ that, in 
this case all correlation functions up to 6 -point which 
are related through Ua ( 1 ) symmetry will be degener¬ 
ate, making the anomalous breaking of Ua ( 1 ) invisible in 
these correlation functions. Thus, if Ua( 1) breaking is fi¬ 
nite through the nondegeneracy of the 2 -point correlation 
functions such as Xn — Xs, the eigenvalue density must 
be nonanalytic in to 2 . Two such possible forms of the 
infrared eigenvalue spectrum, compatible with (ipi/s) = 0 
but x-k — Xs 7 ^ 0 for to —^ 0, have been speculated in 
0 , namely p(X,m) ~ m 2 5(\) and p(X,m) ~ |m|. The 
functional form of the infrared eigenvalue density in the 
chirally symmetric phase of QCD remains an open and 
interesting theoretical question. 

The global Ua( 1) breaking at T = 0 is intimately con¬ 
nected to the presence of topologically nontrivial configu¬ 
rations of the QCD gauge fields | 8 }. It is well known that 
localized topological structures like instantons give rise 
to zero modes of the Dirac operator and the correspond¬ 
ing wavefunctions remain localized 0 - The occurrence 
of near-zero modes can possibly also be traced back to 
the underlying topology of the gauge field configurations. 
For example, the particular form of the eigenvalue den¬ 
sity of the Dirac operator p(X,m) ~ m 2 5(X) in the in¬ 
frared can be motivated from the fact that a small shift 
from zero of the near-zero modes resulting from the weak 
interactions among widely separated instantons and an- 
tiinstantons can be neglected, leading to a <5(A) behavior. 
The to 2 factor naturally arises from the two light fermion 
determinants. At high enough temperatures, it has been 
shown that the dilute gas of instantons is a reasonable 
description of pure gauge theory P . 0 . Within this 
approximation the instanton density is suppressed with 
decreasing value of the gauge coupling and eventually 
vanishes when T —> oo [l 8 |. This dilute gas model is ex¬ 
pected to be a good description of the high temperature 
phase of QCD as well 0 . The ^(l) breaking can be ex¬ 
plained within such a model only for small enough values 
of the gauge coupling, i.e. at sufficiently high tempera¬ 
ture and for small sizes of the instantons. It is however 
unclear whether such a mechanism can explain the Ua(X) 
breaking for the more relevant temperature range from T c 
to a few times T c . Near T c the instantons and antiinstan- 
tons may not be widely separated and weakly interac ting , 
as described by the Instanton Liquid Model (ILM) [201 ]. 
In this model, chiral symmetry breaking arises due to 
the fermion modes associated with strongly interacting 
and overlapping instantons. As the temperature is in¬ 
creased, it was proposed that there is a transition from a 
liquid phase of disordered instantons and antiinstantons 


to a phase of instanton-antiinstanton molecules [U, 0] . 
The chiral symmetry restoration at finite temperature 
may not necessarily be due to the suppression of the in¬ 
stantons (23| but rather due to the temperature depen¬ 
dence of the fermion determinant, which favors polarized 
instanton-antiinstanton molecules. Thus, it is entirely 
possible that for T < T c < 2T C other nonperturbative 
mechanisms responsible for 17,4(1) breaking may gener¬ 
ate an accumulation of near-zero modes leading to a more 
complex form of the infrared eigenvalue spectrum. 

Since topological structures in QCD are inherently 
nonperturbative, lattice QCD techniques are ideally 
suited to address issues related to the Ua{X)- Anoma¬ 
lous 174 ( 1 ) breaking at high temperature was studied on 
the lattice by looking at the nondegeneracy of the 2 point 
correlation functions x-k ~ Xs using staggered fermion for¬ 
mulation, which preserves a remnant of the continuum 
chiral symmetry on the lattice, in [0 and more recently 
with an improved staggered fermion formulation in [25[ . 
In both cases, it was observed that Ua{ 1) was not effec¬ 
tively restored at T > T c . Recently, this issue was re¬ 
visited in a study of the infrared eigenvalue spectrum of 
highly improved staggered quarks (HISQ) 0 and a sim¬ 
ilar conclusion was reached [23]. However, for staggered 
fermions the connection between topology and fermion 
zero modes and hence the reproduction of the index the¬ 
orem is a very subtle issue 0. Attempts to study this 
problem with fermions with exact chiral symmetry on 
the lattice have produced .juxta pos ing results. Studies 
with domain wall fermions [l6l I29f but with a heavier 
than physical pion mass of 200 MeV support the sce¬ 
nario that 17 a( 1) remains broken for T >T C . The eigen¬ 
value density of the Dirac operator for T < 1.2 T c has a 
small peak structure in the infrared favoring the form of 
p(A,m) ~ to 2 <5(A), which can largely account for the ori¬ 
gin of the axial anomaly. However, even in these studies a 
clear separation of the zero and near-zero modes was not 
possible due to moderate residual chiral symmetry break¬ 
ing effects induced by the mixing of the left and right 
handed fermions along the finite fifth dimension. On 
the other hand, another independent preliminary study 
with so-called optimal domain wall fermions and physical 
pion mass reports the restoration of Ua (1) above T c [30:]. 
A study using overlap fermions restricted to the trivial 
topological sector of QCD on relatively small volumes 
also suggests that 1 / 4 ( 1 ) is effectively restored at T ~ T c 
[3l| . However, it is well known that simulations with 
fixed topology are more sensitive to finite volume effects 
and, at present, it is difficult to perform calculations with 
larger lattice volumes or for fluctuating topology in the 
case of overlap fermions due to prohibitively large com¬ 
putational costs. 

Measuring the underlying topology of the gauge fields 
on the lattice requires careful analysis. The gauge fields 
are defined as links connecting the adjacent lattice sites, 
which can be continuously deformed to unity. To study 
localized topological structures one has to remove the 
ultraviolet fluctuations of the fields or approach succes- 
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sively to the minimum of the classical action. The latter 
is done by cooling the gauge configurations (32j. The ul¬ 
traviolet fluctuations can be reduced using smearing [H], 
which involves replacing each gauge link by an average 
over the neighboring links. It is then possible to mea¬ 
sure the topological charge on the lattice using the dis¬ 
cretized version of the integrated FF operator. The most 
popularly used smearin g te chnique, known as the hyper- 
cubic (HYP) smearing (34[, has been shown to provide 
a good estimate of the topological susceptibility. How¬ 
ever, successive smearing may lead to small instantons 
being undetected and a change in the large scale struc¬ 
ture of the gauge fields. Alternatively, one can make use 
of the index theorem (35|, which relates the difference be¬ 
tween the number of right and left-handed fermion zero 
modes to the topological charge of the gauge fields. The 
advantage of this method is that it naturally connects 
topological structures and the eigenmodes of the Dirac 
operator, which are global quantities and depend on the 
gauge links of the entire lattice. 

In the present work we address the temperature depen¬ 
dence of Ua (1) and probable microscopic mechanisms re¬ 
sponsible for its breaking in the high temperature phase 
of QCD by studying the infrared eigenmodes of over¬ 
lap fermions [HI on the background of dynamical (2+1) 
flavors of HISQ gauge field configurations with nearly 
physical fermion masses and large volumes. This HISQ 
discretization scheme has been used for extensive studies 
on QCD thermodynamics jl’, 36| and has small discretiza¬ 
tion errors, resulting in the least taste symmetry breaking 
among all commonly used staggered fermion discretiza¬ 
tions. Preliminary studies with HISQ fermions [371 also 
provide hints that in the continuum limit for two vanish¬ 
ingly small light quark masses, the QCD chiral transition 
may belong to the 3-dimensional 0(4) universality class. 
The issues about the lack of an index theorem of the 
HISQ is overcome by using the Overlap Dirac fermions 
to probe the topology of the HISQ gauge configurations. 
Overlap Dirac fermions circumvent the Nielsen-Ninomiya 
No-go theorem [39[ by sacrificing the ultra-locality crite¬ 
rion, preserve an exact chiral symmetry [381 ] and an exact 


index theorem [401 ] and reproduce the correct anomaly 
41[ even at nonzero lattice spacing. Employing the index 
theorem for the overlap fermions, the topological struc¬ 
tures in SU(2) [jj|] as well as SU(3) pure gauge theories 
[Hill have been studied earlier. However, the rela¬ 
tionship between chiral and Ua( 1) symmetry also cannot 
be addressed within the framework of pure gauge theory 
and the presence of light dynamical fermions is necessary 
to address this question. Furthermore, it is a priori not 
evident whether the same dilute instanton gas picture 
also applies for QCD with near-physical, light dynami¬ 
cal fermions as the presence of light fermions would lead 
to interactions between instantons and induce anomaly 
effects. 

This work is structured as follows: In Sec. El we pro¬ 
vide all necessary computational details pertaining to 
this work. In IIII Al we check the distribution of the topo¬ 


logical charge measured by using the exact index theorem 
of overlap fermions. In IIII Bl we present and discuss our 
results on the eigenvalue distribution of overlap fermions 
on the dynamical HISQ configurations. The contribu¬ 
tion of the low-lying eigenmodes towards Ua{ 1) break¬ 
ing is discussed in IIII Cl while iHTDl contains our results 
about the functional form of the eigenvalue density and 
its implications for Ua{ 1) breaking. In IIII El we verify the 
robustness of the occurrence of near-zero modes. From 
IIII FI to IIII Cl we discuss various different properties of 
the zero and near-zero eigenmodes. Finally, in IIVI we 
summarize and conclude this work. Preliminary results 
of this work were previously presented in [43j. 

II. COMPUTATIONAL DETAILS 

The set of (2 + l)-flavor HISQ configurations used in 
this work was generated by the HotQCD collaboration 
Q. Two lattice sizes were used in this study, 24 3 x 6 and 
32 3 x 8. The strange quark mass m s is set to its phys¬ 
ical value and the light quark mass in all these sets of 
configurations are chosen to be mi = m s / 20, which cor¬ 
responds to a Goldstone pion mass of = 160 MeV in 
the continuum. We studied 5 sets of configurations, two 
at T ~ T c , two at T ~ 1.2 T c and one at T ~ 1.5 T c . Here, 
T c = 154(9) MeV [|| is the chiral crossover temperature 
in the continuum limit. Near T c , in addition, we studied 
configurations generated by the Bielefeld-BNL collabora¬ 
tion [33 with lattice size 32 3 x 6 and a light quark mass 
of mi = m s / 40, which corresponds to m^ = 110 MeV. 
This was to study whether the 1 / 4 ( 1 ) breaking survives 
as the chiral limit is approached. We considered 90-160 
configurations of each set, typically separated by 100 tra¬ 
jectories, and computed the eigenvalues of the overlap 
Dirac operator on them. The lattice sizes, strange to 
light quark mass ratio, temperatures and relevant statis¬ 
tics are shown in Tab. Q] 


A 3 x N t 

mi/rris 

T [MeV] 

N 

N x 

24 3 x 6 

1/20 

162.3 

120 

200 

32 3 x 6 

1/40 

162.3 

90 

400 

OO 

X 

CO 

CM 

CO 

1/20 

165.6 

120 

200 

24 3 x 6 

1/20 

199.0 

100 

100 

CO 

to 

CO 

X 

OO 

1/20 

196.0 

100 

100 

32 3 x 8 

1/20 

237.1 

160 

50 


Tab. I: Lattice size (IV 3 x N T ), mass ratio (mi/m 3 ), tem¬ 
perature ( T ), number of configurations ( N ) and number of 
eigenvalues that were computed per configuration (N\) for 
each ensemble. 

We probe the low-lying eigenmodes of these HISQ 
gauge ensembles through the use of the massless over¬ 
lap Dirac fermion operator 

D ov =M[ 1 + 7 5 sgn [75 D w (-M)]\ , (3) 
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where Dw is the standard Wilson Dirac operator with 
the parameter 0 < M <2. 

For the implementation of the sign function in the 
overlap operator, we computed the lowest 20 eigenvec¬ 
tors of D' w Dw using the Kalkreuter-Simma (KS) Ritz 
algorithm [4J|. The sign function was computed for 
these low modes explicitly, while for the higher modes 
it was approximated by a Zolotarev rational function. 
The number of terms in the Zolotarev function was kept 
to be 15. The overlap operator satisfies the Ginsparg- 
Wilson(GW) relation with a deviation of no more than 
10 ”' at low temperatures and 10 -10 at high tempera¬ 
tures. The square of the sign function deviated from 
identity by about 10 ” 7 - 10 -9 . 

For each temperature, 50 lowest eigenvalues of Dj v D ov 
were computed using the KS algorithm. The zero modes 
of DlD ov come with chiralities ±1. The nonzero eigen¬ 
values come in degenerate pairs with chiralities having 
opposite signs but equal magnitudes, which is usually 
different from unity. These features of the spectrum al¬ 
low us to distinguish between the near and exact zero 
modes within a few iterations of the KS algorithm. The 
KS algorithm was run until the relative error on the 
nonzero eigenvalues of Dl v D ov were estimated to be lower 
than 10” 4 on average and the separation between zero 
and nonzero modes was clearly seen. In most cases, the 
number of eigenvalues was later increased by computing 
the eigenvalues of PD ov P , where P is the projection to 
righthanded or lefthanded modes. This projected opera¬ 
tor has the advantage that it only takes half the time to 
be applied and each nonzero eigenvector can be related to 
a pair of eigenvectors of Dl v D ov , further reducing com¬ 
putation time and memory requirement. However, the 
algorithm becomes quite unstable if the subspace that 
P projects onto contains zero modes and can only be 
used on the opposite chiralities after identifying the zero 
modes. 

For some applications, the eigenvectors of D ov were 
required and not only those of Dl v D ov . Each degener¬ 
ate pair of nonzero eigenvectors of the squared opera¬ 
tor spans a two-dimensional space that also contains two 
eigenvectors of D ov , which are related to each other by 
an application of 75 and have eigenvalues that are com¬ 
plex conjugates of each other. These eigenvectors could 
be obtained by applying an appropriate unitary transfor¬ 
mation to each of the original pairs. 

We also checked the optimal value of the parameter M 
used in the construction of the overlap operator. From 
partially quenched studies it is known that for certain 
choices of M, the corresponding D^Dw can have very 
small eigenvalues, leading to the presence of spurious zero 
modes in the overlap operator [45|. We verified that for 
configurations without zero modes the choice of M did 
not affect the eigenvalues significantly within our preci¬ 
sion. Moreover, for configurations with zero modes we 
chose M such that the sign function and the GW rela¬ 
tion were determined with highest accuracy, ensuring the 
best implementation of the overlap operator. Except for 


a few cases, especially near T Cl we chose M = 1.8. 

For the 32 3 x 6 lattice with a light quark mass of 
mi = m s j 40, the gauge configurations were rough and 
the convergence of the KS for d{j,% was slow, leading 
to imprecise estimates of the GW relation and the sign 
function. I 11 this case, we did two levels of HYP smear¬ 
ing to smoothen out the ultraviolet fluctuations. This 
enabled us to achieve a more precise estimation of the 
overlap sign function, similar to the precision achieved 
for the other ensembles. The effects of the smearing are 
further discussed in Section UlIEI 


III. RESULTS 

A. Topological charge distributions 

The topological charge Q was measured by counting 
the number of zero-modes of the overlap operator and 
determining their chiralities, 

Q = n + - n- , (4) 

where n + (n_) is the number of zero modes with chiral¬ 
ity +1 (—1). Since the underlying HISQ gauge configu¬ 
rations were generated at zero strong CP violating angle 
9 , at any temperature all the topological sectors should 
be spanned. However, the configurations may have been 
trapped in one topological sector and the autocorrelation 
times in such cases may be large. To avoid autocorrela¬ 
tion effects we usually chose the configurations to be sep¬ 
arated by 100 Rational Hybrid Monte-Carlo trajectories. 
The time histories of the topological charge are shown in 
Fig. Q] It is evident that the autocorrelation effects are 
under control. The distribution is ergodic enough and 
on average (Q) ~ 0, for all temperatures. This gives us 
confidence that the statistics is sufficient in our present 
study. Moreover, we also observe many configurations 
with \Q\ > 1 in all studied ensembles. In fact, even at 
1.5 T c more than a third of the total number of configu¬ 
rations have |< 2 | = 1 , confirming the importance of the 
< 2^0 configurations. 


B. Eigenvalue spectra 

In this section we show the eigenvalue density of the 
overlap operator at three temperatures - near T c , at 
1.2 T c and at a yet higher temperature of 1.5 T c . The 
overlap fermion matrix in Eq. <[3j) is a normal matrix. In 
the complex plane, its dimensionless eigenvalues, A, lie 
on a circle centered at M and with a radius M, obeying 
|A — M\ 2 = M 2 . The eigenvalues measured in our study 
lie on the circle very close to the origin with a very small 
real part. Hence, we always plot the eigenvalue density 
as a function of A, where aA = ImA. On the lattice, the 
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Fig. 1: Time histories of the topological charge, calculated 
from the zero modes of the overlap Dirac operator, for the 
HISQ configurations at T ~ T c , T ~ 1.2 T c and T ~ 1.5 T c . 
Configurations belonging to the same production stream are 
connected with lines. 


eigenvalue density is defined as 




q 3 pO) = ^3^ X! ^ ( aA ~ aA 0 > ( 5 ) 

where the sum only includes values on the left part of 
the semicircle with Re A < M, and excludes values near 
A = 2 M even though their imaginary part would also be 
small. 

The eigenvalue distribution at three different temper¬ 
atures is shown in Figs. [2J [3] and [I] The spectrum is 
truncated at some large eigenvalue since we measure only 
a finite number of them. We indicate the point beyond 
which the spectrum is not trustworthy anymore by a ver¬ 
tical line in red in each of the plots. It is estimated by 
first taking the highest computed eigenvalue of each con¬ 
figuration and then taking the minimum of these values 
over all the analyzed configurations. 

As we emphasized earlier, the KS algorithm allowed 
us to distinguish the zero modes from the near-zero 
modes using the chirality properties of the corresponding 
eigenvectors. In general, the eigenvalue distribution has 


Fig. 2: The eigenvalue density of the overlap operator on 
32 3 x 8, 24 3 x 6 and 32 3 x 6 HISQ configurations near T c . 
Only nonzero modes are included. The vertical line denotes 
the range of validity due to the finite number of computed 
eigenvalues. 


three distinct features - the zero mode peak, a near-zero 
mode accumulation and the bulk eigenvalue region. Near 
T c , the first bin contains a large contribution from zero 
modes which are omitted in Fig.[2]to focus on the infrared 
physics of only the near-zero eigenvalues. At this tem¬ 
perature, we do not observe any gap in the infrared part 
of the eigenvalue spectrum. The near-zero modes and 
the bulk modes appear to overlap significantly and the 
near-zero modes tend to develop a peak towards the in¬ 
frared region. This peak becomes sharper as the light sea 
quark mass is lowered from mi = m s /20 to mi = m s /40 
at fixed lattice spacing 1/6T. It also becomes sharper 
when we go to a finer lattice, from N r = 6 to N T = 8 
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X/T 



X/T 


Fig. 3: The eigenvalue density for HISQ configurations using 
the overlap operator at 1.2 T c . The lattice sizes are 32 3 x 8 
and 24 3 x 6, respectively. The red line marks the range of 
validity. 



X/T 



Fig. 4: The eigenvalue density for 32 3 x 8 HISQ configurations 
using the overlap operator for 1.5 T c for all values of Q and 
also separately for the Q = 0 sector. The red line marks the 
range of validity. 


at a fixed pion mass of 160 MeV. This trend suggests 
that the near-zero mode accumulations will remain as 
the chiral and the continuum limits are approached. 

At temperatures 1.2 T c and 1.5 T c , both the zero modes 
denoted by the red bar and the near-zero and bulk modes 
are shown in Fig. [3] and Fig. |4] The separation between 
the near-zero mode accumulation and the bulk eigenvalue 
region becomes even more evident with increasing tem¬ 
perature. At 1.2 T c , we study the eigenvalue spectrum 
at two different lattice spacings to estimate whether the 
infrared part of the spectrum is strongly affected by the 
lattice cutoff effects at higher temperatures. Keeping the 
physical bin size the same in units of X/T for compari¬ 
son, we observe that the infrared region of the eigenvalue 
density remains practically unchanged when the lattice 
spacing goes from 1 /6 T to 1/8 T at a fixed temperature T. 
This gives us confidence that the near-zero modes are not 
due to dislocations of the gauge fields. A more detailed 
study about the lattice artifacts is given in Sec. IIIIDI 

The number of zero and near-zero modes both de¬ 
crease as the temperature is increased to 1.5 T c as shown 
in Fig. []] There is a small peak of near-zero modes, 
while the number of bulk eigenvalues starts to rise very 
slowly and only gives a significant contribution beyond 
Ao — 0.4 T. This is reminiscent of some kind of band 
edge separating the two different regimes of eigenvalues, 


which is studied in detail in Sec. IIII Gl Even at this tem¬ 
perature we do not observe a gap in the infrared sector of 
the eigenvalue spectrum. The presence of these near-zero 
modes is not due to the fact that we are also sampling 
configurations belonging to nonzero topological sectors in 
our study. This is evident from the lower panel of Fig. [4j 
where we show the eigenvalue distribution of only those 
configurations with topological charge Q = 0 at 1.5 T c . 
The presence of near-zero modes is also observed for this 
particular subset of configurations. 


C. Near-zero eigenmodes and axial symmetry 
breaking 


As introduced in Sec. 1 u = Xn ~ X5 as defined in 
Eq. f2j is a measure that quantifies £/4(1) breaking. 
Through a chiral Ward identity it can also be obtained 


w = 


m 


Xconn 


( 6 ) 


from the chiral condensate (ifiip) = y (tr (D rr ^d m D m )'/ 
and the connected chiral susceptibility Xconn = 
y (<9 m tr (D^d m D m )), where D m = D ov (l-am/2M) + 
am is the Dirac operator for overlap quarks with a (va¬ 
lence) quark mass m. 
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Thus, in terms of the eigenvalues of the overlap oper¬ 
ator 


a 2 uj = 


<M> 


N 3 N T [(am) 2 

2(om) 2 (4M 2 — |A| 2 ) 2 


( 7 ) 


(Ett 2 „ , ■ 

\a#o |jA| 2 (4M 2 — (am) 2 ) + 4 (am) 2 M 2 


The first term is the contribution from the zero modes 
which vanishes in the thermodynamic limit. 

Thus, having determined the low-lying eigenvalues of 
the overlap operator, uj can be computed from them. 
However, to explore the physics of the underlying HISQ 
configurations the overlap valence quark mass that en¬ 
ters Eq. 0 has to be tuned against some physical quan¬ 
tity measured on the same gauge configurations. In the 
present work we adopt a very simple strategy: we roughly 
tune the strange valence quark mass, m s , and then study 
w as a function of the light valence quark mass within a 
range of mi = m s / 20 to mi = m s /2. To tune the strange 
valence quark mass we use the renormalized difference of 
the pseudo-scalar ( rj sS ) and the (connected) scalar sus¬ 
ceptibilities in the strange quark sector, m 2 w s /T 4 . We 
calculated this renormalized quantity using strange va¬ 
lence overlap fermions and matched it with the same 
quantity calculated independently for the strange HISQ 
sea quarks. For the overlap fermions this quantity was 
estimated in two parts. We first calculated this quan¬ 
tity from the predetermined low-lying eigenvalues using 
Eq. 0, and then added the contribution of higher eigen¬ 
values by performing inversions on the eigenspace orthog¬ 
onal to the low-lying eigenmodes using random source 
vectors. The result of the strange valence quark mass 
tuning near T c is shown in Fig. [5] In order to monitor 
finite volume effects, we considered the exact zero modes 
separately and do not observe a significant contribution 
from them to this quantity. The tuned masses obtained 
with and without the zero modes differ by about 5%. 

In Fig. 1 we show a renormalized measure of Ua( 1) 
breaking, namely mjm s w/T 4 , for a range of the light va¬ 
lence quark mass between m s / 20 and m s /2. Assuming 
a Breit-Wigner distribution for the near-zero mode peak 
to model the S( A) like distribution discussed in the intro¬ 
duction, i.e. p(X)/T 3 = p 0 A/(A 2 + X 2 ), the contributions 
of the near-zero modes in this renormalized measure of 
I/a( 1) breaking can be characterized as 



Fig. 5: The tuning of the strange valence quark mass for over¬ 
lap fermions on N T = 6 (blue) and N T = 8 (red) lattices near 
T c . The horizontal lines mark the results for m 2 s u s /T A (cf. 
Eq. 0) independently obtained for the HISQ sea fermions. 
Filled (empty) points denote the overlap result with (without) 
the zero mode contribution. 


discussed in Sec. IIII Dl will also give a finite value and 
should be approximated in a partially quenched study. 
For light valence quark masses near to or smaller than 
the smallest near-zero eigenvalues, the computed quan¬ 
tity will approach zero, which is visible for the lowest 
masses near T c . This is a finite volume effect in the sense 
that a larger volume will sample more eigenvalues in the 
near-zero mode region, which reduces the magnitude of 
the smallest eigenvalue and pushes this effect towards 
zero. Beyond that, we observe a smooth dependence on 
mi/m s that is compatible with the Breit-Wigner ansatz 
and which is independent of the lattice spacing. It is also 
evident from Fig. [G] that the contribution of the near¬ 
zero modes to w is substantially larger than that from 
the bulk modes. When going from T ~ T c to T ~ 1.2 T c , 
this quantity does not decrease significantly, supporting 
our conclusion that Ua( 1) is not effectively restored si¬ 
multaneously with the chiral symmetry. 


D. The functional form of the eigenspectra 

In order to understand the general functional form of 
the eigenvalue density and in particular the near-zero re¬ 
gion, we make a fit ansatz consisting of a Breit-Wigner 
peak for the near-zero modes and a polynomial behavior 
for the bulk part of the spectrum of the form, 


mim s A + 2 mi 

—U oc pomsj-r- -rj . 

T 4 (A + mi) 2 


Pi A) 

T 3 


PoA 

A 2 + A 2 


+ cA“ . 


( 9 ) 


In our partially quenched setup, where only the valence 
light quark mass mi is varied, this quantity has a finite 
value in the chiral limit mi —> 0. In a full dynamical 
setup, its behavior will be governed by the dependence 
of A and p 0 on the light sea quark mass. The simplest 
case of both A and po being proportional to the light 
sea quark mass, which is compatible with the trends 


We address three issues in this section. Firstly, a gen¬ 
eral idea about the dependence of the near-zero mode 
peak on the sea quark mass is necessary to understand 
what happens in the chiral limit. Secondly, it is impor¬ 
tant to check the dependence of the near-zero modes on 
the lattice cutoff to establish that these are physical and 
not mere lattice artifacts. Finally, the leading exponent 
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Fig. 6: A renormalized measure of Ua{ 1) breaking for a range 
of valence light quark masses, m s /20 < mi < m s /2 for en¬ 
sembles with different N T at T ~ T c and T ~ 1.2 T c . The 
hlled points denote the contribution from near-zero modes 
(A < Ao), while the empty points were calculated only from 
the bulk modes (A > Ao). 


that characterizes the rise of the bulk also provides in¬ 
formation regarding the restoration of 1/4(1), hence its 
dependence on the temperature and lattice cutoff needs 
to be studied. 

The fit to the eigenvalue spectrum near T c for differ¬ 
ent sea quark masses is shown in Fig. [7] The error bars 
for each bin have been determined by a jackknife proce¬ 
dure over the set of gauge configurations. The param¬ 
eter A , characterizing the width of the near-zero mode 
peak, falls from 0.35(4) T to 0.151(7) T when going from 
mi = m s /20 to m s /40 and the prefactor po, indirectly 
controlling the height of the peak, also goes down from 
0.28(3) T to 0.191(7) T. This generally supports the pic¬ 
ture of a delta function like peak forming in the chiral 
limit, with a decreasing peak height. However, as far as 
one can tell from two data points, the dependence of po 
does not seem to be quadratic in the light sea quark mass 
near T c , arguing against the dilute instanton gas picture 
as a good description of QCD at this temperature. 

To examine cutoff effects one needs to compare a renor¬ 
malized version of the eigenvalue density at different lat¬ 
tice spacings. One way to renormalize the eigenvalues 
is to scale them by the previously tuned strange quark 
mass. The corresponding renormalized eigenvalue den¬ 
sity is then m s p(A) since it leaves the quantity m s (ipijj) 



Fig. 7: Eigenvalue distribution at T ~ T c for N r = 6 and two 
different sea quark masses compared to fits with Eq. ©. 


unchanged under renormalization. We therefore take the 
same ansatz as in Eq. (©) with A —> X/m s and the density 
replaced by its renormalized definition m s p(A)/T 4 . Fits 
to the renormalized spectrum in dimensionless units at 
two different temperatures, T c and 1.2 T c , are shown in 
Fig. El Taking a closer look at the near-zero mode peak, 
it is evident that the accumulation of near-zero modes is 
almost independent of the lattice spacing and is unlikely 
to be just a lattice artifact. 

Finally, we discuss the exponent a of the characteristic 
A“ rise of the bulk eigenvalues. As mentioned in Sec. Q] 
under the assumption of analyticity of the eigenvalue den¬ 
sity in ?n 2 , detailed analytical calculations based on up 
to 4-point chiral Ward identities show 0 that in the 
chiral symmetric phase of QCD the leading A depen¬ 
dence should be similar to that for the free theory, i.e. 
lim m _>o p(A,m) ~ A 3 . In such a case, the effect of Ua( 1) 
breaking should be invisible in at least up to 6-point cor¬ 
relation functions. In light of this, it is interesting to 
characterize the rise of the bulk eigenvalues. As shown 
in Fig. [71 near T c the rise of the bulk eigenvalues for our 
two N t = 6 lattices with quark masses mi = m s /20 and 
mi = m s /40 is described by the exponents a = 0.92(5) 
and a = 0.98(4), respectively. A fit to the renormal¬ 
ized eigenvalue spectrum near T c further yields 0.86(2) 
for N t = 8, see Fig. [8] Thus, near T c a linear rise of bulk 
eigenvalues is favored for both lattice spacings and quark 
masses. Interestingly, chiral perturbation theory [461 ] and 
the ILM [47] show that a contribution to p( A) linear in 
A, which results in a non vanishing connected susceptibil¬ 
ity Xconn, is absent for two light flavors and only present 
for Nf > 2. On the other hand, for staggered fermions 
away from the continuum limit, taste violations lead to 
Xconn ^ 0 also for Nf = 2 HU . 

At 1.2 T c the rise of the bulk eigenvalues has a different 
exponent. From Fig. [S] it is evident that a ~ 2, indepen¬ 
dent of the lattice spacing. Note that a bulk eigenvalue 
density rising quadratically with A does not contribute to 
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Fig. 8: The renormalized eigenvalue spectra near T c (top) and 
at 1.2 T c (bottom), both for two different lattice spacings. The 
lines indicate fits to Eq. © with A — > \/m s and p(X)/T 3 —>• 
m s p(A)/T 4 . 


u>. A similar linear behavior for T ~ T c and a quadratic 
rise for T ~ 1.2 T c of the bulk eigenvalues was also ob¬ 
served in the previous study with domain wall fermions 
|2?)1 with a heavier pion mass, corroborating that the bulk 
rise is independent of the sea quark mass. On the other 
hand, the characteristic free theory like cubic rise of the 
bulk eigenvalue density was only observed at 1.5 T c as 
shown in Fig. [9] At 1.5 T c , the near-zero mode peak re¬ 
duces significantly and the separation between the bulk 
and the near-zero modes is distinctly visible. Our fit 
ansatz for the bulk is modified accordingly as (A — Ao)“ 
to represent this feature, which gives a smaller y 2 per de¬ 
grees of freedom than the original ansatz in Eq. ©. The 
parameter a and the goodness of fit at different temper¬ 
atures are compiled in Tab. QT1 


E. Robustness of the zero and near-zero modes 

Detecting topological objects with fermion zero modes 
has the advantage that by construction the zero modes 
depend on all the gauge links distributed on a lat¬ 



Fig. 9: Eigenvalue distribution at 1.5 T c together with 
a (A — Ao)“ fit to the bulk eigenmodes. 


N 3 x N t 

mi/m s 

T [MeV] 

a 

y 2 /dof 

24 3 x 6 

1/20 

162.3 

0.92(5) 

1.32 

32 3 x 6 

1/40 

162.3 

0.98(4) 

1.84 

CO 

to 

05 

X 

00 

1/20 

165.6 

0.86(2) 

0.92 

24 3 x 6 

1/20 

199.0 

1.9(2) 

1.16 

CO 

to 

05 

X 

oo 

1/20 

196.0 

1.9(1) 

1.21 

CO 

to 

05 

X 

00 

1/20 

237.1 

3.0(4) 

1.30 


Tab. II: Lattice size ( N 3 x N T ), mass ratio ( mi/m s ), temper¬ 
ature (T), the exponent a characterizing the A“ rise of the 
bulk eigenvalues A and the goodness of the fits performed on 
the eigenvalue distribution. 

tice. Still, if the underlying gauge fields are not smooth 
enough, the method might be hampered by the pres¬ 
ence of unphysical fermion modes localized on structures 
called dislocations, which typically have a smaller clas¬ 
sical action than instantons. These are lattice artifacts, 
i.e. effects of finite lattice spacing, and should disappear 
as the continuum limit is approached. It is therefore im¬ 
portant to make sure that the observed infrared fermion 
modes are physical and do not solely arise as lattice ar¬ 
tifacts. If the zero and near-zero eigenvalues are entirely 
due to the presence of dislocations, they are expected to 
disappear as the gauge fields are smoothed using smear¬ 
ing techniques. To check this we performed HYP smear¬ 
ing on the 32 3 x 8 configurations with mi/m s = 1/20 at 
1.5 T c . 

Such smoothing methods are mandatory if one wants 
to compute the topological charge by means of a dis¬ 
cretized version of its field theoretic definition 

Q = sh j A ‘ xF ^ ■ < 10 > 

The optimal number of smearing levels is usually chosen 
such that the topological charge measured this way on the 
smeared configurations has a value close to an integer. 
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Q 


Fig. 10: The distribution of the topological charge, Q, at 
1.5 T c for 32 3 x 8 HISQ configurations measured from the 
zero modes of the overlap operator and also from the gluonic 
operator FF on the same gauge configurations after 10 levels 
of HYP smearing. 


On our present lattices it turned out that 10 levels of 
HYP smearing were sufficient to give an integer value of 
the FF operator summed over the whole lattice at 1.5 T c . 
On all of these smoothed configurations we could verify 
that the value of Q obtained from FF matched exactly 
with the topological charge obtained by counting the zero 
modes of the overlap operator after smearing. 

In Fig. [TOl we compare the histogram for the topologi¬ 
cal charge measured by counting the fermion zero modes 
of unsmeared configurations with that measured using 
the purely gluonic observable FF on the same config¬ 
urations after smearing. Apparently, some zero modes 
disappear in the course of the smearing process; yet, the 
comparison indicates that the fermionic zero modes for 
these configurations do not arise only due to gauge field 
dislocations and are reflecting continuum physics. 

A similar behavior is observed for the near-zero modes. 
In Fig. [lTJ the comparison of the eigenvalue spectrum of 
the overlap operator on the original unsmeared config¬ 
urations with the smeared ones reveals that the num¬ 
ber of near-zero modes is somewhat reduced by smearing 
but that the near-zero mode accumulation is still present 
even after a substantial amount of smearing, indicating 
that these are not mere lattice artifacts like dislocations. 
We also found that the typical eigenvectors associated 
with the near-zero modes on a smeared configuration ap¬ 
pear to be slightly less localized compared to the un¬ 
smeared case, suggesting that the reduction of the near¬ 
zero modes may be caused by the loss of small instantons 
due to smearing. 

Finally, as already shown in Fig. [51 near T c as well as 
at 1.2 T c the comparison of the renormalized eigenvalue 
spectra at two different lattice spacings, 1/6T and 1/8T, 
indicates that the near-zero mode accumulation remains 
nearly unchanged as the lattice spacing is reduced. Thus, 
it is unlikely that the near-zero mode accumulation arises 
primarily due to lattice artifacts. 



Fig. 11: The eigenvalue density at 1.5 T c and N T = 8 before 
and after 10 steps of HYP smearing. The empty boxes are ob¬ 
tained when including zero modes and the dashed lines mark 
the ranges of validity. 


F. Profiles of the zero and near-zero modes at high 
temperature 

The fermion zero mode associated with an instanton is, 
at T = 0, localized in the region occupied by the instan¬ 
ton. At nonzero temperature however, the compactifica- 
tion of Euclidean time leads to periodic copies of instan¬ 
tons. Such classical finite action solutions of the gauge 
fields on the manifold R 3 x S 1 are known as calorons. 
When the instanton size is much smaller than 1/T, the 
copies do not feel the effect of the neighbors and behave 
like zero temperature instantons. However, the overlap 
between the instanton copies may be larger and their 
sizes become comparable to 1/T. Explicit solutions are 
known for trivial [49| as well as nontrivial p5(l [BH ] holon- 
omy. In the case of trivial holonomy it has been observed 
that the caloron turns into a magnetic monopole [52] in 
pure gauge theory. Calorons with nontrivial holonomy 
have more interesting features with monopole substruc¬ 
tures IyI . However, in our studies we have not attempted 
to take a detailed look at the monopoles. 

In order to gain more inside into the structure of the 
infrared modes in QCD, we looked at the profiles of the 
zero and the near-zero modes at 1.5 T c . Examples rep¬ 
resenting the majority of configurations with topological 
charge |<5| = 1 are shown in Figs. IT2l and IT3l In these 
figures the density of the wavefunctions, ip' (x)ip(x), is 
measured along two spacetime directions, summing over 
the other two directions and the internal degrees of free¬ 
dom (color and spin). We observe that the zero modes are 
localized along the spatial as well as the temporal direc¬ 
tions. In a few cases, the width of the zero modes in the 
compact temporal direction is somewhat larger and the 
density profile represents the overlap between the nearest 
copies. 

The near-zero modes typically exhibit a two-peak 
structure, in the density as well as in the chirality 
(ip 1 * (x)j 5 ip(x)) profile, cf. Fig. [131 In the latter case, the 
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Fig. 12: Space-time profile of a zero mode at 1.5 T c for a 
typical gauge configuration with Q = 1. 


chirality contained in the two peaks is of opposite sign. 
These profiles provide a strong hint toward a picture 
where two zero modes with equal and opposite chirali¬ 
ties interact weakly, becoming a pair of near-zero modes 
when properly superposed. 

The fermion zero mode il>n(x) associated with an in- 
stanton is known analytically [171 ]. giving a density of 
the form 

V4(a#o(z) = 2 , l P 2V? , (H) 

+ p z ) J 

where p is the radius of the instanton. When three of the 
spacetime coordinates are integrated over, the density 
along the remaining fourth coordinate, say y, becomes 

= 2 ( 2/2 + ^ 2)372 ' ( 12 ) 

An estimate of the instanton size p can be obtained by 
either finding the distance where this integrated den¬ 
sity falls below l/\/8 of its maximal value or by fit¬ 
ting Eq. (fl2l) to the measured density. When estimat¬ 
ing the size in the temporal direction, the second ap¬ 
proach has the advantage that it can also accommodate 
cases where the size of the instanton ~ 1/T and the pe¬ 
riodic copies in the temporal direction have noticeable 
overlap. This is achieved by replacing the fit function by 
EL - n fr(r + k/T). 


Fig. 13: Space-time profile of a near-zero mode at 1.5 T c for 
the same gauge configuration depicted in Fig. 1121 


At 1.5 T c , we measured the wavefunction density of the 
zero modes along each coordinate direction for all config¬ 
urations with | Q | = 1 by summing over the other three 
directions. The radii of the profiles along the x , y and 
2 directions, p x , p y and p z , were averaged over to give 
a spatial radius of p a = 0.223(8) fm, essentially indepen¬ 
dent of which method was used. Along the temporal di¬ 
rection, using the adjusted fit ansatz we obtained a radius 
that was only slightly larger, namely p T = 0.24(1) fm. 

In order to study the distribution of the distance be¬ 
tween the instanton and antiinstanton forming a pair, 
the lattice points with the lowest and highest chiral den¬ 
sity of the corresponding near-zero modes were identi¬ 
fied and the distance between them was measured. The 
result is shown in Fig. [15] We compare this with the 
expected distribution of the separation if an instanton 
and an antiinstanton were to be distributed on the lat¬ 
tice randomly, independent of each other. It can be seen 
that the measured distribution compares quite well with 
the random distribution with a slight excess at small sep¬ 
arations, further supporting the picture that instantons 
and antiinstantons are weakly interacting. 

In a weakly interacting random ensemble of such topo¬ 
logical objects, the assumption of independent occur¬ 
rence results in a Poisson distribution for the total num¬ 
ber n of instantons and antiinstantons , 

P K (n) = e~ K K n /n \, 


(13) 
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Fig. 14: Space-time profile of chirality of the same near-zero 
mode depicted in Fig. [13] at 1.5 T c . 


0.5fm lfm 1.5fm 2fm 2.5fm 3fm 



Fig. 15: Distribution of distances between the instanton and 
antiinstanton that couple to give rise to a near-zero mode. 
The blue points show the expected distribution if the instan¬ 
ton and antiinstanton were to be distributed randomly and 
independently of each other on a 32 3 x 8 lattice. 


where k is a parameter that is equal to the ensemble av¬ 
erage (n) as well as the variance a 2 and n is obtained 
by counting the number of eigenvalues below a cutoff 
Ao- To fix the cutoff, we fitted Eq. (THU) to the distri¬ 
bution of n for different cutoffs and compared the y 2 
per degree of freedom. The value that was closest to 1, 
namely 1.03, was obtained for Ao = 0.44 T. The result¬ 
ing distribution is shown in Fig. [TB] together with the 



0 2 4 6 8 10 


Fig. 16: Configuration distribution of the total number n of 
zero and near-zero modes at 1.5 T c and Poisson fit to the data. 


Poisson fit which gives k = 4.5(2). At this cutoff, we 
obtained (n) = 4.50(14) and a 2 = 4.2(4) from averaging 
over the configurations, confirming that instantons and 
antiinstantons indeed occur almost independently. 

Taking 4.5(2) as the average total number of instan¬ 
tons and antiinstantons per configuration, we obtain a 
density of 0.147(7) fm~ 4 , which is much lower than the 
value predicted from the strongly interacting dense ILM, 
1 fm~ 4 [23 ■ However, since we are at a temperature much 
higher than the chiral crossover temperature T c , it is not 
surprising that we get values expected rather from a di¬ 
lute and weakly interacting ensemble. 

Altogether, the findings described in this section are 
giving support to a picture in which for sufficiently high 
temperatures of T > 1.5 T c , the infrared behavior of QCD 
can be described as that of a dilute gas of instantons and 
antiinstantons which are interacting weakly with each 
other. 


G. The localization properties of the eigenmodes 


The profiles of the near-zero modes have suggested that 
these are localized structures. To further quantify their 
localization properties and those of the bulk modes, we 
study the so-called participation ratio (PR), defined for 
a normalized eigenvector tp(x) of the Dirac operator as 


PR = 


1 

Nj*N r 





(14) 


It is the fraction of the total lattice volume occupied by 
the eigenmode. If the eigenvector is distributed equally 
on the entire four-volume, this quantity is unity. 

First we use the PR to corroborate our observation 
made in Sec. IIII FI If indeed the near-zero modes repre¬ 
sent weakly interacting instanton-antiinstanton pairs, the 
PR of a typical near-zero mode should be about twice as 
large as that of a zero mode. The comparison of the PR 
of the near-zero modes and of the zero modes of |Q| = 1 
configurations is shown in Fig. 1171 The ratio of the aver¬ 
age PR of a near-zero and that of a zero mode indeed is 
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Fig. 17: Participation ratio of zero modes of \Q\ — 1 config¬ 
urations and near-zero modes (A /T < 0.4) as well as their 
average values as function of the eigenvalue A at T ~ 1.5 T c . 


1.85. The PR values of near-zero inodes fluctuate about 
the mean value, so not all of them support this picture, 
but there is a significant fraction that does. 

At 1.2 T c , the PR histograms are measured for two val¬ 
ues of the lattice spacing and compiled in Fig. [T5] The 
low-lying eigenvalues are more localized than the bulk 
modes; however, the PR gradually increases as one goes 
towards the bulk. At 1.5 T c , the low-lying modes below 
A < 0.4 T average to a PR value of about 0.01. As one 
enters the bulk eigenvalue region for A > 0.4 T, see for 
example Fig.[9l there appears to be a rise of the PR. The 
value of A ps OAT may thus be considered as a mobility 
edge separating the localized near-zero eigenstates from 
the delocalized bulk states. The presence of localized as 
well as delocalized states is observed in disordered semi¬ 
conductors whose dynamics is described by the Anderson 
Hamiltonian. In the Anderson model, the electron states 
at the band edge are localized whereas the states at the 
band center remain delocalized within the lattice. The 
corresponding eigenvalues of the Anderson Hamiltonian 
change from Poisson statistics at the band edge to ran¬ 
dom matrix theory statistics at the band center. 

A comparison with random matrix model predictions 
can be achieved by looking at the distribution of the dis¬ 
tance s between two consecutive eigenvalues, i.e. the level 
spacing distribution. In order to understand its univer¬ 
sal properties, it is necessary to map the eigenvalues onto 
new values using an unfolding procedure I54| , thereby re¬ 
moving the non-universal global scale of the system. By 
construction, the unfolded eigenvalues that are obtained 
by this method have a mean level spacing of unity. 

While the unfolded level spacings of the localized 
near-zero modes should follow a Poisson distribution be¬ 
cause of their mostly independent occurrence, the bulk 
modes should be strongly mixed. If the bulk is highly 
disordered, the corresponding level spacing distribution 
should follow the same distribution as the eigenvalues 
of a random matrix of an appropriate symmetry group. 
The Dirac operator for QCD with matter fields in the 
fundamental representation falls into the same symme- 



A/T 


Fig. 18: The PR for eigenvectors at 1.2 and 1.5 T c . The 1.2 T c 
data are for two different lattice spacings to study the cutoff 
dependence of our results. 


try group as Gaussian Unitary ensembles (GUE) whereas 
for the case of two colors it is in the same universal¬ 
ity class as Gaussian Orthogonal ensembles (GOE). The 
level spacing distribution for the unfolded bulk eigenval¬ 
ues of the overlap operator on HISQ configurations at 
1.5 T c is shown in Fig. 1191 It shows an agreement with 
a random matrix theory with Gaussian unitary matri¬ 
ces. This is in general agreement with similar studies 
of the localization of low-l ying modes done before for 
the quenched theory M [56H58II and also with dynamical 
staggered fermions HIM f6oj ~ on smaller lattice sizes. 


IV. CONCLUSIONS 

In this work we have investigated the temperature de¬ 
pendence of the anomalous Ua( 1) symmetry breaking 
in the high temperature phase of QCD with two light 
quark flavors. To this end we have employed the over¬ 
lap Dirac operator exploiting its property of preserving 
the index theorem even at nonvanishing lattice spacing. 
We have applied the overlap operator on large volume 
HISQ gauge field configurations and computed its low- 
lying eigenmodes. We observe the presence of zero as 
well as near-zero modes in the investigated temperature 
range of T c < T < 1.5 T c . By comparing the low-lying 
eigenmodes from two lattice spacings and studying the 
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Fig. 19: The universal level spacing distribution for A > 0.4 T 
at 1.5 T c compared to the Random Matrix theory Gaussian 
Unitary ensemble (GUE) and the Gaussian Orthogonal en¬ 
semble (GOE). 


effects of smearing we have shown that these infrared 
modes are not mere lattice cutoff effects. 

We mainly analyzed configurations which have been 
obtained at a light sea quark mass corresponding to a 
pion mass of 160 MeV. However, within the set of con¬ 
figurations at our disposal, at a temperature near T c we 
could confirm the accumulation of the near zero eigenval¬ 
ues also at a quark mass considerably below its physical 
value. 

By quantifying the contribution of the near-zero eigen- 
modes to a specific combination of two point correla¬ 
tion functions, ~ Xs- we conclude that these modes 
are primarily responsible for the anomalous breaking 
of the axial symmetry in QCD still being visible for 
T c T < 1.5 T c . Through detailed studies of their 


spacetime profiles, localization properties and distribu¬ 
tions over gauge configurations we have shown that for 
T ~ 1.5 T c the near-zero modes follow the behavior as 
expected of a gas of widely separated, weakly interacting 
instantons and antiinstantons. At 1.5 T c we find the den¬ 
sity of (anti)instantons to be 0.147(7) fm -4 , with a typi¬ 
cal radius of 0.223(8) fm. At this temperature, the spatial 
volume of our lattice was ~ (3.3 fm) 3 with 1/T ~ 0.83 fm, 
suggesting that the instanton gas is indeed dilute, the 
instanton size is smaller than 1/T and our chosen vol¬ 
ume being large enough to accommodate more that one 
instanton-antiinstanton pair on average. In conclusion, 
our study suggests that at T ~ 1.5 T c the origin of global 
Ua(X) breaking in QCD is due to the dilute gas of weakly 
interacting instantons and antiinstantons. 

For an independent confirmation of our results, it 
would clearly be desirable to carry out a similar anal¬ 
ysis with dynamical chiral fermions. While the lattice 
spacing effects in this work appear to be small, it will 
further be necessary to control the subtle extrapolations 
to the continuum as well as the chiral limit in such a 
future investigation. 
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